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1.  Project  Summary 


The  goal  of  this  project  is  to  improve  our  scientific  understanding  of  dispersion  issues 
over  complex  terrain  and  publish  results  in  scientific  journals.  Three  scientific  articles 
and  one  student’s  Master  thesis  are  under  preparation  based  on  the  results  supported  by 
this  project. 

The  first  paper,  entitled  “The  role  of  advection  of  fluxes  on  modeling  dispersion  in 
convective  boundary  layers”,  is  aimed  at  examining  the  importance  of  the  advection 
terms  for  dispersion  in  a  convective  boundary  layer.  An  Eulerian  three-dimensional 
higher-order  closure  dispersion  model  is  presented.  The  model  uses  mean  wind  and 
turbulence  values  from  a  second  order  atmospheric  boundary  layer  model.  The  dispersion 
model  is  validated  against  results  from  tank  and  field  experiments  and  compared  to 
results  from  Lagrangian  dispersion  models.  The  results  show  good  agreement  with 
experiment  and  Lagrangian  modeling  results  for  point  source  dispersion  in  a  convective 
boundary  layer.  Sensitivity  studies  of  the  model  helped  to  identify  the  roles  that  advection 
and  horizontal  transport  terms  in  the  equations  for  the  fluxes  play  in  simulating  the 
essential  features  of  pollutant  dispersion.  The  results  from  the  sensitivity  tests  show  that 
the  famous  features  of  dispersion  from  a  point  source  in  the  convective  boundary  layer  - 
with  an  ascending  plume  during  ground  level  release  and  descending  for  a  lifted  point 
source  -  is  caused  by  the  advection  term  in  the  equation  for  the  vertical  flux.  Furthermore, 
it  is  shown  that  there  is  a  tendency  for  the  plume  to  split  horizontally,  which  is  also 
caused  by  the  advection  term  in  horizontal  fluxes. 

The  article  entitled  “Modeling  of  atmospheric  pollutant  dispersion  in  complex  terrain  - 
Development,  validation  and  comparison  of  a  higher-order  closure  mesoscale  model”  is  a 
two  part  paper.  In  part  I,  data  from  the  extensive  field  program  MOHAVE 
(Measurements  Of  Haze  And  Visibility  Experiment)  was  used  to  evaluate  the 
meteorological  output  of  two  atmospheric  mesoscale  models,  the  higher-order  MIUU 
model  (Meteorological  Institute  of  Uppsala  University’s  model)  and  the  MM5  (the  fifth- 
generation  Penn  State/NCCAR  Mesoscale  Model).  Part  II  presented  a  higher-order 
closure  dispersion  model  suitable  for  modeling  dispersion  in  complex  terrain.  The  results 
of  this  model  have  been  compared  against  tank  experiment  as  well  as  the  MOHAVE  data. 

In  the  student’s  Master  thesis  a  trajectory  model  is  developed.  This  model  has  two- 
dimensional  (on  an  isentropic  surface)  and  three-dimensional  options.  This  model  has 
many  applications  such  as  tracing  the  dust  trajectory  from  the  Gobi  desert. 

The  content  of  the  first  paper  will  be  presented  in  Section  2.  The  other  two  papers,  which 
are  still  under  preparation,  will  be  reported  in  Section  3.  The  thesis  will  be  reported  in 
Section  4. 

2.  The  role  of  advection  of  fluxes  on  modeling  dispersion  in  convective 
boundary  layers 
2.1  Introduction 


In  recent  years,  due  to  an  increased  interest  in  environmental  problems,  considerable 
attention  has  been  focussed  on  means  of  predicting  concentration  levels  in  turbulent 
boundary  layers.  In  practice,  predictions  are  often  required  for  dispersion  in  complicated 
flows  over  coastal  areas  or  other  complex  terrain,  since  air  pollution  sources  are  often 
situated  in  such  areas.  To  obtain  the  flow  fields  over  complex  terrain  one  usually  solves  a 
complete  set  of  hydrodynamic  and  thermodynamic  equations.  This  provides  the 
meteorological  parameters  needed  for  the  dispersion  calculations.  The  dispersion 
calculations  follow  different  lines; 

One  type  of  approach  is  to  use  the  standard  Gaussian  dispersion  model,  modified  to 
account  for  terrain  influences.  Examples  of  such  models  are  described  by  Egan  (1984), 
Hunt  et  al.  (1979),  Mistra  (1980),  Andren  (1987),  and  Enger  (1990c).  With  this  approach, 
terrain  adjustment  factors  are  introduced  to  correct  for  the  calculation  of  standard 
deviation  of  the  Gaussian  concentration  distribution.  The  weakness  of  this  approach  is 
that  the  adjustment  factors  are  site-dependent.  Furthermore,  most  of  these  types  of 
models  are  only  applicable  on  the  idealized  topography  on  which  they  are  formulated. 


Another  approach  is  to  use  Lagrangian  particle  dispersion  models.  Atmospheric 
dispersion  in  a  Lagrangian  model  is  simulated  by  tracing  a  large  set  of  particles  driven  by 
wind  and  turbulence  fields  predicted  by  the  atmospheric  model.  Subsequent  positions  of 
each  particle,  representing  a  discrete  element  of  pollutant  mass,  are  computed  from  the 
relations.  The  resolvable  scale  components  of  wind  velocity  are  obtained  directly  from 
the  mesoscale  meteorological  model,  and  turbulent  components,  are  generated  with  the 
aid  of  schemes  such  as  a  Markov  chain  scheme  or  a  fully  random  walk  scheme.  Examples 
of  such  models  are  described  in  Uliasz  (1994)  and  Lamb  (1979). 

A  third  approach  is  to  use  an  Eulerian  diffusion  model,  which  starts  from  the  mass 
continuity  equation.  A  first-order  closure  model  for  complex  terrain  was  used  by  Segal  et 
al.  (1982)  to  study  dispersion  in  the  Greater  Chesapeake  Bay  area.  Enger  (1986)  has 
shown  the  capability  of  a  second-order  closure  dispersion  model  in  simulating  observed 
dispersion  features  in  a  convective  boundary  layer.  The  model  was  also  used  to  simulate 
dispersion  in  complex  terrain  in  Enger  (1990b)  and  Enger  and  Koracin  (1995).  The 
present  model  is  a  further  development  of  this  model.  The  present  model  solves 
prognostically  the  mean  concentration,  C,  and  the  second  order  moments,  which  include 
concentration,  i.e.  ^  and  .  In  the  two-dimensional  version  of  the  model 

used  in  Enger  (1986),  the  equations  for  the  horizontal  fluxes  (m'c' and  v'c')  were 
neglected.  In  Enger  (1990b)  a  three-dimensional  version  of  the  model  was  used,  that 
version  all  tendency  equations  for  the  second-order  moments  were  included,  but  with  the 
horizontal  advection  terms  and  the  horizontal  transport  terms  neglected.  Close  to  the 
source,  these  terms  have  been  shown  to  be  as  important  as  other  terms  in  the  equations 
(Fackrdl  and  Robins,  1982).  In  Enger  (1990b)  these  terms  were  neglected  by  supposing 
the  horizontal  advection  and  horizontal  dispersion  balanced  each  other.  The  aim  of  the 
present  study  is  therefore  to  examine  the  importance  of  these  terms  for  dispersion  in  a 
convective  boundary  layer. 


The  basic  equations  and  parameterizations  used  in  the  dispersion  model  are  presented  in 
section  two  and  three.  The  experiments  are  presented  in  section  four.  The  model 
validation  is  discussed  in  section  five  and  the  sensitivity  tests  in  section  six.  Conclusions 
from  the  present  study  are  given  in  section  seven. 

2.2.  The  basic  model  equations 

The  present  model  is  an  Eulerian  diffusion  model,  which  starts  from  the  mass  continuity 
equation  that  in  tensor  notation  reads 


dt  ^  dXj  dXj 


(1) 


For  explanation  of  symbols  see  nomenclature.  The  corresponding  equations  for  the 
turbulent  fluxes,  if  we  neglect  molecular  terms  and  the  effect  of  the  Coriolis  force  on  the 
covariances,  are 


rhi'  r'  chi'  c'  - - dU  du'M'.c'  1  dv'  J?  - 

Sll±  =  -U  .--u'  c'^ - - Lc'2E--1j-0<c'  (2) 

dt  '  ax,  '  ^  8x.  dx.  p,  dx,  ©0 


The  equation  for  the  covariance,  c'G',  is 


dc'G'  rr  dc'G'  -r^ac  -r-;a0  du'jC'G' 

- =  -C7  . - u\0 - M  ,  c' - - - 

dt  ^  dx,.  dx,.  dx.  dx. 


D 


(3) 


The  pressure  covariance  terms  in  equation  (2)  are  parameterized  according  to  Enger 
(1986)  as 


1  -c'-^  =  a,  c'--^c'e\ 


Pq  ax,  X 


(4) 


where  ai  is  a  constant,  and  A,  is  a  master  turbulent  length  scale. 

In  the  model  version  used  in  Enger  (1986)  -  a  2-dimensional  version  -  the  horizontal 
transport  terms  were  neglected,  and  in  Enger  (1990b)  the  advection  terms  were  also 
neglected  in  the  equations  for  the  horizontal  fluxes.  According  to  measurements  by 
Fackrell  and  Robins  (1982)  these  terms  are  as  important  as  other  terms,  close  to  the 
source.  Therefore,  in  the  present  model  the  terms  are  retained  and  the  consequence  of 
neglecting  any  of  the  terms  shall  be  discussed  in  the  sensitivity  tests. 

The  vertical  transport  term  is  parameterized  by  use  of  a  gradient  diffusion  approximation, 
according  to  Donaldson  (1973)  and  Mellor  (1973),  as 


(5a) 


dM>'u/c'  ^  d 

and  the  horizontal  transport  term  is  parameterized,  in  the  same  manner,  as 


dx^  dx^  dXf 


i=l,2 


(5b) 


where  pz  and  pH  are  constants,  which  shall  be  determined  later.  Moreover  we  treat  the 
transport  terms  in  the  equation  for  the  covariance  c'O'  in  a  similar  fashion,  i.e. 


dc'e' 


^  =  .  i-1.2 

dz  dz  dz  dXj  dx^  St,. 


(6) 


The  molecular  destruction  term  D  is  parameterized  according  to  Lumley  (1975)  as 


D  =  a,^c'e' 

X 


(7) 


where  a2  is  a  constant. 

The  values  for  the  constants,  ai,  aj,  Pz  and  Ph  are  0.3465,  0.144,  0.36  and  0.36, 
respectively.  The  formulation  of  the  master  length  scale,  X,  is  given  in  the  Appendix  B. 
The  atmospheric  variables  are  obtained  from  the  higher-order  closure  mesoscale  model 
developed  at  the  Meteorological  Institute  of  Uppsala  University,  the  MIUU  mesoscale 
model.  The  full  description  of  the  model  is  not  repeated  here,  but  the  main  three- 
dimensional  equations  are  listed  in  the  Appendix  A. 

2.3  Numerical  method  and  boundary  conditions 

A  finite-difference  numerical  method  is  used  to  solve  the  set  of  dispersion  equations.  The 
model  uses  an  expanding  telescoping  grid  mesh  with  its  origin  at  the  source.  This  means 
that  we  get  a  denser  grid  mesh  near  the  source  where  the  plume  is  narrower  than  further 
downstream.  In  all  experiments  in  the  study  we  used  71x71  grid  points  in  the  honzontal 
and  31  grid  points  in  the  vertical.  A  grid  mesh  distance  of  50m  x  50m  x  10m  at  the  source 
and  a  model  dimension  of  50  km  x  10  km  x  5  km  are  used.  The  prognostic  equations  in 
the  dispersion  model  are  solved  by  using  an  advection  scheme  that  is  of  third-order  both 
in  space  and  time  (Enger  and  Grisogono,  1997).  The  diffusion  is  solved  with  a  semi- 
implicit  scheme  with  weight  0.75  on  the  future  time  step.  The  Arakawa  staggered  C-gnd 
is  used  for  the  model  structure.  The  following  boundary  conditions  are  applied  at  the  top 
of  the  model 


dC  dw'c' 


:«'c'  =  v'c'  =  c'^'  =  0 


dz  dz 


(8) 


and  the  following  are  applied  at  the  surface 


dC  -r-,  ^  ^  ^  ^ 

- =  w'c  =  — - —  =  — —  =  — 

dz  dz  dz  dz 


(9) 


2.4  Experiments 

Various  experiments  were  designed  with  the  aim  at  validating  and  investigating  the 
sensitivity  of  the  dispersion  model.  The  one-dimensional  version  of  the  MIUU  model  was 
used  to  simulate  an  ideal  convective  boundary  layer  (CBL).  Figure  1  shows  the  profiles 
of  some  of  the  atmospheric  parameters  in  the  simulated  CBL.  These  atmospheric 
variables  were  spread  into  the  3-dimensional  dispersion  model  to  obtain  a  horizontally 
homogeneous  CBL  in  the  model.  The  boundary  layer  depth  obtained,  as  the  height  with 
the  minimum  heat  fiux,  is  about  846  m,  see  Figure  Id.  Since  we  intend  to  compare  the 
model  results  with  those  of  the  Deardorff  and  Willis  (1975),  and  Willis  and  Deardorff 
(1978,  1981)  tank  experiments,  we  need,  according  to  Deardorff  and  Willis  (1975),  to 
ensure  that  the  simulated  CBL  satisfies  the  criteria  1.5^*  <  U  <6w,.  From  Figure  1  we 
can  see  that  the  mean  wind  is  about  2.3  m  in  the  CBL.  The  convective  velocity  scale  is 
1.1m  s'^  hence  the  criteria  is  satisfied.  With  this  atmospheric  condition  a  CONTROL 
experiment  was  earned  out  Avith  the  full  version  of  the  model  for  validations.  A  series  of 
other  experiments  were  also  carried  out  to  investigate  the  model  sensitivity  to  some  of  the 
terms  in  equation  (2)  and  (3).  The  model  set-up  for  each  of  these  sensitivity  tests  will  be 
discussed  later. 

2.5  Model  validation 

The  CONTROL  experiment  was  designed  to  validate  the  dispersion  model  against  results 
from  tank  experiments  for  surface  and  elevated  sources.  Therefore  we  will  simulate 
dispersion  from  point  sources  at  0.062,  0.25,  0.5  and  0.75  of  the  boundary  layer  height. 
For  ease  of  comparison  with  laboratory  and  field  experiments,  the  results  are  presented  in 
dimensionless  form.  The  dimensionless  concentration  of  pollutants  and  the  dimensionless 
crosswind-integrated  concentration  are  defined  as 

and  =  T  ,respectively  (10) 

Q  y  Q 

The  horizontal  y-direction  is  normal  to  the  mean  wind  direction  along  x-direction.  The 
dimensionless  downwind  distance,  crosswind  distance  and  height  are  expressed  as 

X  =  ;  Y  =  -^  and  Z  =  —  ,respectively.  (1 1) 

ZiU  Zj  zj 


2.5. 1  Crosswind  integrated  concentration 

We  start  by  looking  at  a  vertical  cross  section  of  the  dimensionless  crosswind  integrated 
concentration,  Cy,  for  surface  and  elevated  sources.  Figure  2a-d. 

For  the  near  ground  level  release  (0.062zi),  shown  in  Figure  2a,  the  plume  centerline  first 
moves  horizontally  with  the  surface  mean  wind  then  rises  rapidly  upward  (at  distance  of 
about  0.4  ZjU/w.)  toward  the  inversion  base  as  it  is  being  advected  by  the  mean  wind. 
The  plume  axis  reaches  an  elevated  height  of  0.9  Zi  at  distance  of  around  2.2  ZjU/w* .  It 
features  a  surface  minimum  around  2.2  ZiU/w, .  The  plume  centerline  also  starts  to 
descend  after  striking  the  inversion  layer  above  to  produce  a  surface  maximum 
concentration  at  distance  4.5  ZiU/  w, .  After  a  distance  of  about  6  ZiU/>r»  the  pollutants  are 
well  mixed  throughout  the  entire  boundary  layer  with  a  dimensionless  concentration  of 
about  0.9,  which  means  that  10%  has  entered  the  stable  boundary  layer  above.  The  main 
behaviors’  are  similar  to  those  in  the  results  from  the  tank  experiment  presented  in 
Deardorff  and  Willis  (1975)  and  the  Eulerian  model  results  in  Enger  (1986).  In  the  tank 
experiment  results  the  plume  centerline  started  to  rise  at  0.5  ZjU/w,  to  an  elevated 
maximum,  which  is  located  roughly  at  the  height  0.8  Zi  and  at  a  distance  of  1.5  ZiU/w, . 
The  surface  minimum  concentration  is  located  around  1.6  ZjU/w*.  In  the  Enger  (1986) 
results,  the  elevated  maximum  is  located  at  height  0.85  Zj  at  a  distance  of  2.5  ZiU/w,  and 
the  minimum  being  located  at  1 .6  ZiU/  w, . 

The  result  for  a  source  at  0.25  z,  is  given  in  Figure  2b.  The  plume  has  a  tendency  to  split 
into  two  parts,  one  rising  and  the  other  one  descending.  The  descending  part  intersects  the 
ground  at  a  distance  of  about  0.5  ZiU/  w. ,  after  which  it  rises  to  merge  with  the  first  part. 
The  plume  then  moves  towards  the  inversion  base.  The  plume  reaches  its  maximum 
height  0.9  Zi  at  a  downwind  distance  of  2.0  ZiU/w, .  There  are  three  maximum 
concentration  features  in  this  case;  two  at  the  surface,  located  at  distances  0.5  ZiU/w,  and 
4.5  ZiU/w.;  the  third  is  an  elevated  maximum  located  below  the  inversion  layer  at  a 
downwind  distance  of  about  2.0  ZjU/ w', .  After  a  distance  of  6  ZjlJ/  w* ,  the  concentration 
is  well  mixed  within  the  entire  CBL  with  concentrations  of  about  0.95.  The  results  are  in 
very  good  agreement  with  those  from  the  tank  experiments  by  Willis  and  Deardorff 
(1978),  the  Lagrangian  model  simulation  by  Lamb  (1978)  and  the  Eulerian  model 
simulation  by  Enger  (1986).  In  their  results,  the  plume  central  line  intercepts  the  ground 
at  0.4,  0.6,  and  0.5  ZiU/w.,  respectively.  The  present  model  gives  almost  the  same  value 
of  both  surface  and  elevated  maximum  concentration  as  observed  in  the  tank  experiment 
results.  It  is  also  important  to  note  that  the  model  simulates  the  ‘rebound’  of  the  plume  at 
the  surface  (i.e.  the  rise  of  the  plume  after  striking  the  ground  surface)  quite  well.  This 
feature,  which  is  well  captured  by  the  present  model  (as  observed  in  the  tank 
experiments),  was  not  simulated  by  the  previous  Eulerian  model  version  (Enger,  1986). 

Figure  2c  shows  a  similar  result,  but  for  a  source  at  the  middle  of  the  CBL.  Like  in  the 
previous  case  with  the  source  at  0.25zi,  the  plume  splits  into  two  parts.  The  descending 
part  intercepts  the  ground  at  about  1.0  ZiU/w. .  The  interception  took  place  at  a  downwind 


Figure  2.  The  dimensionless  crosswind  integrated  concentration  for  point  sources  at  heights 
0.062zi  (a),  0.25zi  (b),  0.5zi  (c)  and  0. 75zi(d)  as  simulated  in  the  CONTROL  run. 


distance  of  about  0.8  ZiU/w.  in  both  the  Willis  and  Deardorff  (1981)  tank  experiment  and 
Lamb’s  (1982)  simulation.  We  also  notice  the  rising  of  the  plume  after  intercepting  the 
ground,  which  is  in  good  agreement  with  the  tank  experiment.  Another  interesting  feature 
is  the  appearance  of  a  minimum  in  source  height  at  a  distance  0.9  ZjU/w,  downstream. 
These  simulated  features,  for  a  source  in  the  middle  of  the  CBL,  are  in  excellent 
agreement  with  the  CONDOR  field  experiment  results  (Moninger  et  al.,  1983)  and  Willis 
and  Deardorff  (1981)  tank  experiment.  The  minimum  normalized  concentration  values 
are  0.75,  0.8,  and  0.6  in  the  present  simulation,  the  tank  experiment,  and  the  field 
experiment,  respectively. 

For  a  source  at  0.75  Zi,  shown  in  Figure  2d,  the  descending  plume  intercepts  the  ground  at 
a  distance  of  about  1.5  ZiU/w,.  This  agrees  well  with  Lamb’s  (1978)  simulation,  but  in 
the  Enger  (1986)  simulation  the  intercept  took  place  at  1.9  ZiU/w.  Also  in  this  case,  the 
present  model  results  show  a  rebound  of  the  plume  after  intercepting  the  ground. 

Therefore  it  can  be  seen  that  the  dimensionless  crosswind  integrated  concentrations  fi-om 
the  present  model  simulations  are  in  excellent  agreement  with  results  fi-om  the  tank 
experiments,  field  experiments  and  Lagrangian  model  simulations.  In  some  extent  the 
present  model  simulations  give  a  better  representation  of  these  features  (i.e.  closer  to 
those  in  the  tank  experiments)  than  the  simulations  of  Enger  (1986). 

Simulated  dimensionless  crosswind  integrated  concentrations,  Cy,  at  ground  level  (for  a 
ground  level  release  at  0.025  Zi)  are  presented  in  Figure  3a.  These  results  are  compared 
with  the  Prairie  Grass  data  (Nieuwstadt,  1980),  the  results  from  a  water-tank  experiment 
(Deardorff  and  Willis,  1975)  and  from  Lamb’s  Lagrangian  model  simulation  ^.amb 
1982).  From  the  figure  it  can  be  noticed  that  the  present  model  results  agree  well  with  the 
above-mentioned  results,  except  in  the  vicinity  of  the  source  where  the  present  model 
diverges  from  the  Prairie  Ghass  data.  This  is  caused  by  the  fact  that  in  the  Prairie  Grass 
data  the  source  is  located  very  close  to  the  surface,  whereas  in  the  present  simulations  the 
source  is  at  21  m  height  (zi  =  846m),  which  results  in  smaller  values  of  surface 
concentrations  close  to  the  source.  There  is  a  minimum  of  the  crosswind  integrated 
concentration  at  around  2  ZiU/w»,.  Afterwards  it  starts  to  increase  and  level  off  to  an 
almost  constant  value  further  downstream  when  the  pollutant  is  well  mixed  throughout 
the  entire  CBL  (Convective  Boundaary  Layer).  The  minimum  corresponds  to  the 
minimum  concentration  at  ground  for  a  source  at  the  height  0.062  Zj,  shown  in  Figure  2a. 

2. 5. 2  Mean  ground-level  concentration 

Now  let  us  examine  the  corresponding  dimensionless  ground-level  concentrations 
CN(x,y,0),  which  are  presented  in  Figure  4.  One  of  the  most  striking  features  in  this  figure 
is  that  the  surface  concentration  distribution  across  the  plume  is  non-Gaussian.  This  non- 
Gaussian  distribution  is  especially  pronounced  for  a  source  at  height  0.062  Zj,  Figure  4a, 
where  a  double  maximum  across  the  plume  is  present.  The  double  maximum  appears 
between  a  downwind  distance  of  about  0.5  and  2.0  ZiU/w* .  This  feature  is  also  present  in 
Lamb’s  (1979)  simulations  with  the  Lagrangian  model.  Such  a  feature  has  also  obseiwed 
in  some  dispersion  measurements  on  the  East  Coast  of  the  US  under  CBL  conditions 
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Figure  4.  The  dimensionless  ground-level  concentration  in  x-y  plane  for  point  sources  at 
heights  ft  062  Zi  (a),  0.25  Zj  (b),  ft.  5  z/  (c)  and  ft.  75  Zj  (d). 


(Nappo  2001,  personal  communication).  As  will  be  discussed  later  in  the  sensitivity  test, 
advection  and  horizontal  diffusion  terms  in  the  equations  for  the  fluxes  play  a  significant 
role  in  the  creation  of  this  feature.  Figure  5  shows  the  ground  level  concentration  across 
the  plume  at  a  distance  of  1.0  ZiU/w.  for  source  heights  0.062zi  and  0.25zi  .  A  similar 
shape,  with  two  maximum  features  in  the  ground  level  concentration,  appears  for  both 
source  heights,  but  is  less  pronounced  for  a  source  at  height  0.25zi. 


2.6  Sensitivity  tests 

The  explanation  for  the  shown  features  of  dispersion  from  a  point  source  in  the  CBL  is 
according  to  Lamb  (1982)  that  the  distributions  of  the  vertical  winds  are  non-Gaussian. 
Instead  there  should  be  short  periods  -  or  small  areas  -  with  high  up-winds  in  the 
convective  cells,  followed  by  long  periods  -  or  large  areas  -  with  slow  descending  air. 
However,  the  Eulerian  model  that  has  been  used  here  does  not  include  any  probability 
density  functions  of  the  vertical  velocity,  but  still  gives  very  similar  dispersion 
distribution  from  a  point  source  in  the  CBL.  The  following  sensitivity  tests  will  show  that 
there  exists  another  explanation  to  the  plume  behaviour  in  a  CBL.  The  sensitivity  tests 
wUl  also  explain  some  horizontal  features  that  can  be  noticed  in  the  simulation  results 
shown  above. 

The  primary  aim  with  the  sensitivity  tests  is  to  investigate  the  relative  importance  of  the 
horizontal  advection  and  diffiision  terms  in  equations  (2)  and  (3).  For  this  purpose  a 
number  of  sensitivity  experiments  were  carried  out  with  model.  In  the  first,  both  the 
advection  and  diffusion  terms  are  neglected.  In  the  second,  only  the  diffusion  terms  are 
neglected,  whereas  in  the  third,  only  the  advection  terms  are  excluded.  Furthermore,  tests 
have  been  performed  to  investigate  the  sensitivity  of  the  coefficient,  Dh,  in  the  diffusion 
terms. 

2. 6. 1  No  Advection,  No  Horizontal  Diffusion  Experiment  (NOAD) 

The  experimental  set-up  for  this  experiment  is  the  same  as  for  the  control  experiment, 
except  that  both  the  advection  and  horizontal  diffusion  terms  in  the  equation  for  all  the 

fluxes  (M^'and  VO')  are  excluded.  It  should  be  noted  that  for  a  steady  state  solution  with 
no  advection  and  no  horizontal  diffusion  terms  in  equation  (2),  the  expression  for  the 
horizontal  crosswind  flux,  v'c',  reduces  to 


-  v'c' 


a,q  5y  dy 


In  our  case,  with  horizontally  homogenous  conditions  the  turbulence  is  constant  at  a 
certain  height,  giving  an  exchange  coefficient  Kh  that  is  constant  for  that  height.  A 
constant  Kh  (Fickian  diffusion)  should  give  a  Gaussian-like  distribution  of  concentrations 
across  the  plume.  In  fact,  the  results  of  the  model  for  this  case  show  that  the  horizontal 
concentration  distribution  is  Gaussian  (see  Figure  6,  circles).  Moreover,  all  the  essential 
dispersion  features  for  a  point  source  in  a  CBL,  i.e.  descending,  ascending  and  ‘rebound’ 
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Figure  6.  Comparison  of  the  lateral  crosswind  profile  of  the  dimensionless  ground-level 
concentration.  Source  height  0. 062zi  and  distance  1. 0  ZjUAv*. 


of  the  plume  and  the  non-Gaussian  distribution  across  the  plume,  disappear.  Figure  7 
shows  the  simulated  dimensionless  crosswind  integrated  concentration  and  the 
dimensionless  ground-level  concentration  as  an  example  of  this  behavior.  Results  from 
this  experiment  show  that  the  advection  and  diffusion  terms  in  the  equation  for  w'c'  are 
crucial  in  simulating  the  vertical  distribution  of  the  concentration  from  a  point  source  in  a 
CBL. 

2.6.2  Only  Advection  Experiment  (OADV) 

As  pointed  out  by  Fackrell  and  Robins  (1982),  the  advection  terms  and  the  diffusion 
terms  are  quite  large  close  to  the  source.  They  also  found  that  the  two  terms  are  almost 
equal  but  with  opposite  sign.  Our  hypothesis  is  that  the  non-Gaussian  distribution  and 
splitting  of  the  plume  is  caused  by  the  fact  that  the  advection  term  in  the  equations  for  the 
fluxes  is  larger  than  the  difibision  terms.  If  we  exclude  the  diffusion  term  but  keep  the 
advection  in  the  equations  for  the  second-order  moments,  the  splitting  of  the  plume 
should  increase  considerably.  Figure  8a  shows  examples  of  this  behavior  for  a  ground 
level  source  (0.062  Zi).  The  results  are  compared  with  the  corresponding  results  obtained 
from  the  CONTROL  experiment.  Figure  4a.  In  Figure  6  the  ground  level  concentration 
values  across  the  plume  at  distance  0.7  ZiU/w,  are  shown  for  the  different  runs:  control 
run  (squares),  NOAD-experiment  (circles),  OADV-experiment  (triangles),  and  a 
simulation  with  higher  diffusion  (constant  Ph  =  1)  than  in  the  CONTROL  run  (stars). 

The  influence  of  the  advection  terms  on  the  dispersion  features  can  be  explained  by 
looking  at  the  fluxes.  If  we  have,  for  example,  a  Gaussian  distribution  (like  in  NOAD 
experiment)  the  horizontal  crosswind  flux  is  zero  at  centerline  and  increases  in  magnitude 
outwards  to  a  maximum  value,  after  which  it  decreases  again.  This  maximum  value  is  of 
course  larger  in  magnitude  closer  to  the  source  than  further  downstream.  The  advection 
term  does  not  effect  the  fluxes  at  the  centerline,  as  they  are  zero  there,  but  away  from  the 
centerline  the  advection  term  causes  these  higher  fluxes  closer  to  the  source  to  move 
downstream.  Which  in  turn  increases  the  flux  value  downstream  and  sharpens  the  flux 
gradient  from  the  centerline  and  out.  Larger  gradients  and  larger  ma^itudes  of  the  fluxes 
cause  an  increased  movement  of  concentration  out  from  the  centerline,  which  causes  the 
plume  to  split.  The  advection  of  fluxes  also  explains  the  noticed  distribution  of 
concentration  in  the  vertical  direction,  with  a  plume  that  split  into  an  ascending  and  a 
descending  part. 


2. 6. 3  Only  Diffusion  Experiment  ( ODIF) 

The  model  set  up  for  this  experiment  is  the  same  as  for  the  NOAD  experiment,  except 
that  the  diffusion  terms  are  included  equation  (2)  and  (3).  The  model  results  replicate 
those  of  the  NOAD  experiment,  where  the  CBL  features  are  absent  and  the  lateral 
crosswind  concentration  profile  is  more  Gaussian-like.  The  diffusion  terms  in  the  second 
order  moment  equations  act  to  smooth  out  the  gradients  of  the  fluxes,  which  means  that 
the  diffusion  terms  act  in  the  opposite  direction  of  the  advection  terms.  If  the  advection 
term  is  larger  than  the  diffusion  term  the  plume  can  split  both  horizontally  and  vertically 
as  shown  in  the  CONTROL  run. 
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Figure  7.  The  model  results  from  NOAD  experiment  for  point  source  at  height  0. 5  Zi.  (a)  the 
dimensionless  crosswind  integrated  concentration  and  (b)  the  dimensionless  ground-level 
concentration  and  for  a  point  source  at  height  0. 062  zt.  (c)  the  dimensionless  crosswind 
integrated  concentration  and  (d)  the  dimensionless  ground-level  concentration  source  at  height 
0.5  Zj. 


Figure  8.  (aj  Dimensionless  ground-level  concentratio 
simulated  in  WOADV  experiment  (b)  dimensionless  grou 
at  height  0. 062  Zj  as  simulated  in  control  run  but  with 


2.6.4  Effect  of  Horizontal  Diffusion  Coefficient,  fin 

This  experiment  was  designed  to  investigate  the  effect  of  the  diffusion  coefficient,  Ph,  on 
the  model  simulations.  The  model  set-up  is  the  same  as  in  the  CONTROL  experiment 
except  that  Ph  =  1,  which  means  a  much  larger  diffusion  of  the  fluxes.  The  results  show 
the  same  features  as  the  CONTROL  experiment  but  without  the  horizontal  splitting  of  the 
plume,  see  figure  8b.  The  constant  Ph  =  1  causes  less  concentration  to  move  out  from  the 
centerline  resulting  in  a  narrower  plume  and  higher  concentration  in  the  center  (see 
Figure  6).  The  model  constant  pH  ==  0.36,  used  in  the  CONTROL  run  was  chosen  to  be 
equal  to  the  vertical  constant  Pz,  as  the  main  turbulence  features  are  caused  by  eddies 
going  through  the  whole  CBL  and  the  time  scales  (and  length  scales)  should  be  the  same 
vertically  and  horizontally. 


2. 7  Summary  and  conclusions 

A  higher  order  closure  dispersion  model  suitable  for  the  airflow  over  complex  terrain  is 
developed,  tested  and  validated  under  convective  conditions.  The  study  has  clearly 
illustrated  that  the  present  higher-order  closure  dispersion  model  is  capable  of  modeling 
dispersion  in  convective  conditions  quite  realistically.  There  have  been  discussions  in  the 
literature  as  to  whether  a  Eulerian  diffusion  model  would  enable  the  prediction  of 
dispersion  in  a  convective  boundary  layer  or  not.  However,  the  results  from  the 
simulations  show  very  good  agreement  with  Willis  and  Deardorff  s  experiments  and  with 
Lamb’s  Lagrangian  model.  The  explanation  for  the  shown  features  of  dispersion  from  a 
point  source  in  the  CBL  is  according  to  Lamb  (1982)  that  the  distributions  of  the  vertical 
winds  are  non-Gaussian  and  have  a  negative  mode  in  the  CBL  and  that  the  probability 
density  profiles  change  with  height.  However,  the  Eulerian  model  that  has  been  used  here 
does  not  include  any  probability  density  fijnctions  of  the  vertical  velocity,  but  still  gives 
very  similar  dispersion  distribution  from  a  point  source  in  the  CBL.  The  sensitivity  tests 
performed  showed  that  there  exists  another  explanation  to  the  plume  behavior  in  a  CBL. 

Sensitivity  studies  with  the  model  helped  to  identify  the  role  the  advection  and  diffusion 
terms  in  the  equations  for  fluxes  play  in  simulating  the  dispersion  from  point  sources  in 
CBLs.  To  simulate  the  horizontal  and  vertical  features  in  the  CBL,  the  advection  terms  in 
the  equations  for  the  second  order  moments  that  involve  concentration  are  shown  to  be 
very  important.  The  horizontal  concentration  as  well  as  the  vertical  concentration 
distribution  are  sensitive  to  the  balance  between  the  advection  and  diffusion  terms  in  the 
second  order  moment  equations. 

The  main  findings  in  the  present  paper  are: 

•  The  descending  and  ascending  of  a  plume  in  a  convective  boundary  layer  is  caused  by 
the  advection  of  the  vertical  concentration  flux,  w'c' . 

•  The  advection  of  the  horizontal  fluxes,  v'c'  and  i/'c',  causes  the  plume  to  split 
horizontally. 


•  The  diffusion  terms  in  the  equations  for  the  fluxes  act  in  the  opposite  direction  of  the 
advection  terms,  but  the  advection  and  diffusion  terms  do  not  balance  each  other 
completely  in  the  CBL. 

3.  Development,  validation  and  comparison  of  a  higher-order  closure  mesoscale 
model 

3.1  Comparison  of meteorological  fields  from  MIUU  and  MM5 

3.1.1  Introduction 

Dispersion  models  need  realistic  three-dimensional  atmospheric  fields  as  input  for 
transport  and  dispersion  studies  of  atmospheric  pollutants.  In  this  study,  data  fi'om  the 
extensive  field  program,  MOHAVE  (Measurements  Of  Haze  And  Visibility  Experiment) 
(e.g.  Green,  1999;  Pitchford  et  al.,  1999),  have  been  used  to  evaluate  the  meteorological 
output  of  two  atmospheric  mesoscale  models:  the  Meteorological  Institute  of  Uppsala 
University’s  model  (MIUU  model)  and  the  Fifth-Generation  Penn  State/NCAR 
Mesoscale  Model  (MM5). 

Project  MOHAVE  was  conducted  in  the  southwestern  United  States  (see  Figure  9)  and 
included  meteorological,  chemical  and  visibility  measurements  at  several  stations.  Upper- 
air  observations  were  made  at  several  sites  using  radiosondes,  wind  profilers  and  a  radio 
acoustic  sounding  system  (RASS).  The  main  goal  of  project  MOHAVE  was  to  assess  the 
effects  of  the  Mohave  Power  Plant  (MOPP),  a  large  coal-fired  power  plant  in  the 
Colorado  River  Valley  in  southern  Nevada,  circa  120  km  to  the  south-southeast  of  Las 
Vegas,  on  visibility  in  the  southwestern  United  States,  in  particular  at  Grand  Canyon 
national  park  (Green,  1999). 

3.1.2  Model  Description 

The  MIUU  model  is  a  higher-order  turbulence  closure  model  of  level  2.5  (according 
to  Mellor  and  Yamada,  1974).  It  applies  a  terrain-following  coordinate  system  (Pielke, 
1984),  is  hydrostatic,  and  includes  prognostic  equations  for  the  mean  variables  and 
turbulent  kinetic  energy  (e.g.,  Enger  et  al.,  1993).  The  turbulence  closure  scheme  was 
corrected  to  include  wall  effects  in  the  redistribution  terms  for  the  second-order  moments 
(Andren,  1990)  in  order  to  account  for  the  influence  of  the  underlying  surface.  Andren 
(1990)  pointed  out  that  this  is  of  importance  particularly  for  air  pollution  applications. 
Additionally,  the  model  includes  routines  for  clouds,  radiation  and  soil  surface 
temperature  and  humidity.  The  model  output,  i.e.,  the  velocity,  temperature  and  turbulent 
kinetic  energy  fields,  is  then  used  in  the  dispersion  model  (see  Section  3.2)  together  with 
the  turbulent  length  scale,  to  compute  all  the  second-order  moments,  that  are  required  in 
the  prognostic  equations  for  the  atmospheric  pollutant  or  tracer  concentrations  (Andren, 
1990). 

A  new  lateral  boundary  condition,  based  on  a  modified  version  of  Davies  (1976)  flow 
relaxation  scheme,  is  used  in  the  model  to  introduce  the  large-scale  synoptic  fields  from 
the  NCEP/NCAR  Reanalysis  project  (Kalnay  et  al.,  1996).  A  terrain-following  flow 


relaxation  parameter  is  applied  in  order  to  allow  mesoscale  phenomena  in  the  planetary 
boundary  layer  to  develop  in  the  vicinity  of  the  lateral  boundaries.  Therefore,  the  flow- 
relaxation  parameter  K  (according  to  Davies,  1976)  is  chosen  to  be  zero  at  the  ground 
surface  and  approaching  its  maximum  value  of  0.0015  at  500  meters  above  ground  level 
(AGL)  for  the  horizontal  wind  components,  and  4000  meters  AGL  for  temperature  and 
humidity.  Between  these  heights  and  the  model  top,  the  parameter  K  is  set  to  0.0015. 
Between  ground  level  and  500  meters  AGL,  or  4000  meters  AGL,  respectively,  a  linear 
interpolation  is  used  for  the  parameter  K.  The  function  F  (according  to  Davies,  1976)  was 
chosen  to  be  zero  everywhere,  and  no  relaxation  was  applied  for  the  vertical  wind 
component.  The  pressure  fields  of  the  NCEP/NCAR  Reanalysis  project  were  interpolated 
to  the  model  grid  points  to  obtain  the  large-scale  synoptic  pressure  gradient  forces. 

Since  Reanalysis  data  were  available  only  every  six  hours,  the  large-scale  synoptic 
pressure,  temperature,  humidity,  and  velocity  fields  were  interpolated  linearly  in  be^een 
these  times.  The  method  allowed  the  model  to  run  for  more  than  9  days,  without 
experiencing  any  instability. 

The  MM5  mesoscale  model  is  well-known  and  thoroughly  described  in  Grell  et  al. 
(1995).  A  detailed  description  is  omitted  here  for  brevity  reasons. 

3.1.3  Model  Simulations 

The  MIUU  model  was  run  fi-om  5  Standard  Mountain  Time  (SMT  =  UTC  -7  hours)  on 
August  5,  1992,  to  5  SMT  on  August  14, 1992.  NCEP/NCAR  Reanalysis  data  have  been 
used  to  prescribe  the  large-scale  synoptic  pressure-fields  and  to  update  the  boundary 
values  of  the  mesoscale  model  using  the  method  described  above.  Surface  elevation  and 
land-use  information  have  been  obtained  from  the  International  Geosphere-Biosphere 
Program  Data  and  Information  System  (Belward,  1996).  Total  cloudiness  was  set  to  0.1 
during  the  entire  period  and  5  K  were  added  to  the  NCEP/NCAR  Reanalysis 
temperatures,  to  get  a  better  agreement  with  measurements.  The  model  was  run  with  30 
vertical  levels,  and  without  using  any  cloud  routines. 

MM5  model  simulations  for  MOHAVE  have  been  performed  and  described  by  Isakov 
(1998).  A  brief  description  of  the  simulations  can  be  found  in  Koracin  et  al.  (2000).  MM5 
was  run  fi-om  17  SMT  on  August  5,  1992,  to  16  SMT  on  August  14,  1992.  MM5 
preprocessing  includes  an  advanced  objective  analysis  of  the  synoptic  data  from  the 
global  network,  providing  detailed  initial  and  boundary  conditions  for  the  simulations.  A 
non-hydrostatic  version  with  35  vertical  levels  was  used.  To  include  more  upper-air 
observations  in  the  initialization  process,  an  expanded  grid  of  60  km  beyond  the  lateral 
boundaries  of  the  model  was  used  for  the  preprocessing. 

In  both  models,  a  horizontal  resolution  of  3  km  was  used.  The  model  domain  covered  the 
area  between  34.5°  and  37°  North  and  112°  to  116°  West  (see  Figure  9).  For  the  MM5 
runs,  terrain  was  obtained  from  a  30”  Defense  Mapping  Agency  Database,  and  then 
considerably  smoothed.  The  terrain  used  in  the  MM5  simulations  does  not  exactly 
correspond  to  Figure  9. 
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Figure  9.  Model  domain  as  used  for  MUU-model  simulations.  The  total  domain  is  used  for  the 
atmospheric  model  simulations,  the  inner  square  is  used  for  the  dispersion-model  simulations. 
Stations  are  as  follows:  MOPP  —  Mohave  Power  Plant  (213  m  ASL);  COCO  -  Cottonwood 
Cove,  West  (274  m  ASL);  COCE  =  Cottonwood  Cove,  East  (201  m  ASL);  DOSP  =  Dolan  Springs 
(853  m  ASL);  DOSl  =  Dolan  Springs,  Tab  (1015  m  ASL);  MEAD  =  Meadview  (905  m  ASL); 
OVER  =  Overton  (424  m  ASL);  OVBE  =  Overton  Beach  (396  m  ASL);  SPMO  =  Spirit  Mountain 
(1498  m  ASL);  SQMO  =  Squaw  Mountain  (1981  m  ASL);  TRUX  =  Truxton  (1350  m  ASL)  and 
LASN  =  Las  Vegas  McCarran  International  Airport  (664  m  ASL).  Circles  indicate 
meteorological  surface  observations,  triangles  meteorological  upper-air  observations  and 
squares  tracer  measurements.  For  details  and  coordinates  see  Pitchford  et  al.  (1999). 


3.1.4  Model  Results 

Results  from  both  models  show  that  the  strong  channeling,  exerted  on  the  flow  by  the 
Colorado  River  valley,  with  heights  of  about  200  -  300  m  ASL  at  the  bottom  of  the 
valley  and  ridges  of  up  to  1800  m  ASL  on  both  sides  of  the  valley,  is  dominating  most  of 
the  time.  Since  the  valley  axis  is  roughly  south  to  north,  the  flow  in  the  valley  is  usually 
to  the  South  in  winter  and  to  the  North  in  summer  (e.g..  Green,  1999).  Figure  10  shows 
an  example  of  the  flow  pattern  at  140  m  ASL  obtained  with  the  MIUU  model. 

Enger  et  al.  (1993)  found,  when  studying  the  flow  regime  of  this  part  of  the  Colorado 
River  valley,  that  the  wind  veering  simulated  over  the  range  of  topographic  elevations  is 
often  larger  than  100  degrees,  in  some  cases  as  large  as  180  degrees.  This  was  also  seen 
in  the  results  of  both  models  and  in  the  measurements. 

3.1.5  Comparison  and  Statistical  Evaluation  against  Surface  and  Upper  Air 
Measurements 

Figure  1 1  shows  a  comparison  of  MIUU  and  MM5  model  results  against  surface 
measurements  at  Mohave  Power  Plant  (MOPP).  In  both  models,  the  agreement  between 
measured  and  modeled  wind  speed  is  quite  poor.  Both  measured  and  simulated  wind 
directions  clearly  show  the  channeling  effect  of  the  Colorado  River  valley.  The  model 
simulated  wind  direction  fields  are  compared  poorly  against  the  measurements  during 
larger  time  periods.  However,  the  MIUU  model  seems  to  perform  slightly  better  for 
temperature,  whereas  for  humidity,  MM5  gives  slightly  better  results. 

Table  1.  Root  mean  square  differences  ofMMS  and  MIUU  model  simulations  compared  to 
upper-air  observations  at  stations  MOPP,  COCE,  DOSl,  MEAD,  OVER  and  TRUX.  Stations  are 
indicated  by  triangles  in  Figure  9.  Number  of  profiles  is  indicated  as  Nprojues-  IF  indicates 
comparison  to  wind-profiler  data,  RS  comparison  to  radio  soundings.  All  differences  are 
calculated  as  the  mean  for  the  entire  height  interval  of  the  profiles,  mostly  from  ground  surface 
up  to  5000  m  AGL. 
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Table  1  shows  a  comparison  of  the  model-simulated  data  with  upper-air  measurements. 
Only  the  time  period  where  results  from  both  models  were  available  was  chosen  for  the 
comparison.  Root  mean  square  (RMS)  differences  (e.g.  Pielke,  1984)  were  calculated 
using  all  available  upper-air  data.  RMS  differences  were  calculated  as  the  mean  over  the 
entire  height  interval  of  the  respective  profile  at  each  station.  Finally,  the  mean  of  the 
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Figure  10.  Horizontal  wind  vectors  at  140  m  AGL  as  simulated  by  the  MIUU  model  at  12  SMT 
on  August  5,  1992.  A  scaling  arrow  of  10  m/s  is  drawn  in  the  lower  left  corner. 


RMS  differences  of  all  available  profiles  at  one  station  was  calculated  and  included  in 
Table  1 .  The  number  of  profiles  used  is  given  as  Nprofiies  in  the  table.  As  can  be  seen  fi'om 
the  table,  both  models  perform  comparatively  equally,  giving  a  RMS  difference  of 
slightly  over  3  m/s  for  the  vector,  around  3  m/s  for  the  u-wind  and  around  4  m/s  for  the  v- 
wind  RMS  differences.  The  comparison  with  RASS  measurements  (indicated  by  asterisks 
in  the  table),  however,  should  not  be  taken  seriously,  since  the  quality  of  these 
measurements  seems  to  be  rather  poor. 

3. 1. 6  Conclusions 

Two  mesoscale  models,  the  MIIJU  model  and  the  MM5  model,  have  been  compared 
to  surface  and  upper-air  observations.  Statistical  methods  give  an  insufficient  answer  to 
the  question,  which  model  actually  is  performing  better.  Thus,  as  Koracin  et  al.  (2000) 
suggested,  tracer  measurements  should  be  used  to  evaluate  wind  fields  from  different 
models.  Studies,  using  both  MIUU  and  MM5  model  results  as  input  to  a  higher-order 
closure  dispersion  model,  will  be  presented  in  Section  3.2. 

3.2  Dispersion  comparisons  between  models 

3.2.1  Introduction 

Many  of  our  air-quality  related  problems  are  concerned  with  pollutant  transport  and 
diffusion  in  complex  terrain.  However,  modeling  pollutant  dispersion  in  such  areas  is  a 
difficult  task  because  of  the  associated  complex  flow  and  turbulence.  The  conventional 
Gaussian  and  Eulerian  K-theory  models  are  inadequate  for  modeling  dispersion  in 
complex  terrain.  Alternatively,  Lagrangian  particle  models  and  higher-order  closure 
models  are  usually  used  for  this  purpose,  e.g.  Yamada  (1985)  and  Enger  and  Koracin 
(1995). 

This  section  reports  a  higher-order  closure  dispersion  (HOCD)  model,  suitable  for 
modeling  dispersion  in  complex  terrain.  A  brief  description  of  the  model  is  presented  in 
3.2.2.  The  model  has  been  applied  to  simulate  dispersion  in  a  convective  boundary  layer 
(CBL).  Comparison  of  the  results  with  those  from  tank  experiments  is  done  in  3.2.3.  The 
model  has  also  been  used  to  simulate  tracer  transport  and  dispersion  in  the  complex 
terrain,  discussed  in  Section  3.1.  The  model  results  are  compared  with  available 
measurements  in  3.2.4. 


3.2.2  Model  Description 

The  present  model  is  an  Eulerian  diffusion  model,  which  starts  from  the  mass  continuity 
equation 


dt  ■'  dXj  dXj 


(1) 


The  corresponding  equations  for  the  turbulent  fluxes,  if  we  neglect  molecular  terms  and 
the  effect  of  Coriolis  forces  on  the  covariance,  are 
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The  equation  for  the  covariance,  c'0' ,  is 

(3) 

di  '  dXj  '  dx^.  dXj.  dXj 

The  pressure  covariance  terms  in  equation  (2)  are  parameterized  according  to  Enger 
(1986)  as 

=  (4) 

p  ^x^  A  3  ^0 

where  ai  is  a  constant,  and  X  is  a  master  turbulent  length  scale. 

The  vertical  transport  term  is  parameterized  using  a  gradient  diffusion  approximation 
according  to  Donaldson  (1973)  and  Mellor  (1973).  The  molecular  destruction  term  D  is 
parameterized  according  to  Lumley  (1975)  as 

D  =  a2^Je' 

A 

The  values  for  the  constants,  ai,  and  aa,  are  0.3465  and  0.144,  respectively.  The 
formulations  of  the  master  length  scale,  X,,  are  described  in  Enger  (1990).  A  finite- 
difference  numerical  method  is  used  to  solve  the  set  of  model  equations.  The  model  uses 
an  expanding  telescoping  grid  mesh  with  its  origin  at  the  source.  This  means  that  we  get  a 
denser  grid  mesh  near  the  source,  where  the  plume  is  narrower  than  further  downstream. 
The  prognostic  equations  in  the  dispersion  model  are  solved  using  a  third-order  advection 
scheme,  described  in  Enger  and  Grisogono  (1997).  The  diffusion  is  solved  with  a  semi- 
implicit  scheme  with  weight  0.75  on  the  future  time  step.  The  Arakawa  staggered  C-grid 
is  used  for  the  model  structure.  The  following  boundary  conditions  are  applied  at  the  top 
of  the  model: 

=  =  =  O  (6) 

dz  dz 

and  the  following  applied  at  the  surface; 
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3.2.3  Model  Simulations  in  Horizontally  Homogenous  CBL  • 

The  model  is  used  to  simulate  the  characteristics  of  plume  dispersion  in  a  CBL  with  point 
sources  at  0.62,  0.25  and  0.5  of  the  boundary  layer  height.  A  one-dimensional  version  ot 
the  MIUU  mesoscale  model  was  used  to  simulate  an  ideal  CBL  used  in  the  dispersion 
model  The  CBL  height,  obtained  as  the  height  of  minimum  heat  flux,  is  about  846  m, 
with  a  mean  wind  speed  of  about  2.3  m/s.  The  convective  velocity  scale  is  calculated  to 
1  1  m/s.  The  CBL  is  similar  to  those  used  in  the  Deardorff  and  Willis  tank  expenments 
(Deardorfif  and  Willis,  1975;  Willis  and  Deardorff,  1978,  1981).  It  satisfies  the  same 

criteria,  i.e.,  1.5w,  <  U  <  6w, . 


Contour  plots  (in  the  x-z  plane)  of  the  dimensionless  crosswind  integrated  concentrating 
simulated  by  the  model  with  different  source  heights  are  presented  in  Figure  12.  T  e 
results  are  compared  with  results  from  the  tank  experiments.  We  can  see  that  the 
simulations  show  the  significant  features  of  dispersion  in  CBL,  as  observed  m  the  t^ 
experiments,  remarkably  well.  The  model  captures  the  rise  and  fall  of  the  plume.  T  e 
‘rebound’  of  the  plume  after  striking  the  surface,  or  the  inversion  layer  above,  is  equally 
well  simulated.  The  maximum  and  minimum  concentration  values,  including  their 
locations,  agree  very  well  with  the  results  from  the  tank  experiments. 

3  2.4  Model  Simulation  in  Complex  Terrain 

The  model  was  then  used  to  simulate  pollutant  transport  and  dispersion  at  Colorado  River 
valley  during  the  period  of  an  intensive  field  program,  the  Measurement  of  Haze  and 
Visibility  Experiment  (MOHAVE).  The  topography  of  the  area  with  the  domain,  used  in 
the  presented  dispersion  simulation,  is  shown  in  Figure  9.  One  of  the  main  objectives  of 
the  field  experiment  was  to  investigate  and  identify  the  possible  short-  and  long-term 
impacts  of  atmospheric  pollutants  from  major  urban  areas  and  industnal  sources 
(especially  from  Mohave  Power  Plant,  MOPP)  on  the  Grand  Canyon  and  its  vicinity.  In 
the  field  experiment,  PTF  ortho-perfluordimethylcyclohexane  (oPDCH)  was  released 
continuously  from  MOPP  during  50  days  in  summer.  Forty-five  percent  of  the  oPDOT 
consists  of  the  isomer  ortho-cis  (oc)PDCH,  which  was  measured  at  the  receptor  sites.  The 
tracers  are  inert,  non-depositing,  and  non-toxic  chemicals.  The  stack  height  of  Ae  power 
plant  is  150  m.  The  effective  source  height  is  calculated  to  be  about  300  m  AGL.  Thus, 
for  the  present  study  an  effective  source  height  of  300  m  AGL  was  used.  The  average 
release  rate  of  oPDCH  for  the  period  is  about  0.042  g/s.  Details  on  the  experiment  and 
measurements  can  be  found  in  Green  (1999)  and  Isakov  (1998).  The  measurements  used 
for  model  comparison  in  the  present  study  are  obtained  from  Isakov  (1998),  but 
converted  from  femtoliters  per  liter  (fl/L)  to  microgram  per  unit  volume  (pg/m  ). 

Dynamic  flow  and  turbulence  fields  over  the  study  area  were  simulated  with  the  MIW 
model  for  the  period  5/8/92  -  14/8/92.  The  results  are  discussed  and  compared  with 
measurements  in  Section  3.1.  For  the  present  study,  only  the  results  of  the  first  three  days 
were  used.  For  use  in  the  dispersion  model,  the  MIUU  results,  which  are  stored  at  1-hour 
intervals,  were  first  interpolated  to  five-minute  intervals.  Thereafter,  the  results  were 
interpolated  to  the  dispersion  model  grid  points.  Figure  13  (a)  and  (b)  show  the 
interpolated  horizontal  wind  vectors  at  100  m  AGL  at  1200  MST  (Mountain  Standard 
Time)  on  5/8/92  and  1700  MST  on  6/8/92.  Note  that  the  wind  vectors  are  plotted  only  at 


Figure  12.  Contours  (in  x-z  plane)  of  the  dimensionless  crosswind  integrated  concentration 
simulated  by  the  HOCD  model  (left)  and  those  from  Deardorff  and  Willis  (right)  for  source 
heights  ofO.  062  zt,  0. 25  z,  and  0. 5  z,.  Sources  are  indicated  by  arrows  on  the  ordinate. 
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Figure  13.  Model  simulations  in  complex  terrain  (Colorado  River  valley)  during  the  MOHAVE 
field  program  in  summer  1992.  Panel  (a)  and  (b)  show  horizontal  wind  vectors  at  100  m  AGL 
interpolated  from  the  MIUU  model  results  at  1200  MST  on  5/8/92  and  0700  MST  on  6/8/92, 
respectively.  Maximum  vector  length  is  13  m/s  and  19  m/s,  respectively.  Panel  (c)  and  (d)  show 
the  HOCD  daily  mean  of  the  simulated  surface  concentration  at  0700 MST  on  6/8/92  and  7/8/92, 
respectively.  Terrain,  receptor  locations,  and  mean  concentration  measurements  are  shown  in  the 
background.  The  values  of  both  simulated  and  measured  concentrations  are  in  KT^  pg/m\ 


every  fifth  grid  point.  Terrain  and  receptor  sites  are  shown  in  the  background.  The 
results  generally  show  southerly  and  southwesterly  flow  in  the  area.  The  southerly  flow  is 
stronger  along  the  Colorado  River  valley  and  follows  the  bending  shape  of  the  valley. 
Effects  of  anabatic  and  katabatic  winds  along  the  mountain  slopes,  due  to  horizontal 
pressure  gradients  created  by  temperature  differences,  are  visible  in  the  simulated  flow  at 
1200  and  0700  MST,  respectively,  especially  around  SPMO  and  the  adjacent  mountains. 

Tracer  release  was  represented  in  the  dispersion  model  by  a  uniformly  distributed  area 
source,  extending  over  the  nearest  four  horizontal  and  three  vertical  grid  points  from  the 
source  location.  The  horizontal  grid  domain  used  is  81x81  km,  with  1  km  resolution  close 
to  the  source.  In  the  vertical,  we  use  3 1  grid  points,  with  50  m  resolution  close  to  the 
release  height.  Simulated  release  started  at  1200  MST  on  the  first  day  and  continued  until 
1200  MST  on  the  third  day.  The  daily  mean  simulated  surface  concentrations  at  0700 
MST  on  the  second  and  third  day  are  presented  in  Figure  13(c)  and  (d),  respectively. 
Terrain,  receptor  sites,  and  mean  measured  concentrations  at  0700  MST  on  the  two  days 
are  shown  in  the  background.  The  values  of  both  simulated  and  measured  concentrations 

are  in  units  10'^  tig/m^. 

From  Figure  13  we  can  see  that  the  simulated  concentrations  are  on  the  same  order  of 
magnitude  as  the  measured  ones.  The  difference  between  the  concentrations  could  be  due 
to  the  fact  that  the  release  of  pollutants  in  the  model  did  not  start  at  the  same  time  as  in 
the  field  experiment.  Also,  insufficient  spatial  resolution  of  both  the  mesoscale  and  the 
dispersion  model  can  affect  the  model  results.  Nevertheless,  the  agreement  between 
simulated  and  measured  concentrations  is  quite  good.  The  model  was  able  to  reproduce 
the  spatial  pattern  of  the  concentrations  over  the  domain.  The  model  results  show  that, 
during  stable  conditions,  the  pollutants  are  transported  aloft  a  few  kilometers  downwind 
before  reaching  the  ground  surface.  During  unstable  conditions,  however,  the  pollutants 
reach  the  ground  surface  in  the  vicinity  of  the  release  point  due  to  stronger  vertical 
mixing.  Moreover,  the  model  simulations  show  that  part  of  the  tracer  released  from 
MOPP  is  transported  to  the  Grand  Canyon  and  its  vicinity. 

3.2.5  Conclusions 

A  higher-order  closure  model  suitable  for  dispersion  in  complex  terrain  is  developed.  The 
model-simulated  features  in  CBL-dispersion  with  different  source  heights  have  been 
discussed.  The  results  are  in  very  good  agreement  with  those  from  tank  experiments. 
Also,  the  model  is  applied  over  a  complex  terrain  near  the  US  Southwest  Coast,  to  study 
transport  and  diffusion  of  tracers  released  at  MOPP  during  the  MOHAVE  1992  field 
experiment.  The  simulated  pattern  and  magnitude  of  concentrations  are  in  good 
agreement  with  measurements.  The  model  results  also  show  that  part  of  the  pollutants 
released  at  MOPP  is  transported  to  Grand  Canyon  and  its  vicinity.  The  simulations  will 
be  refined  in  the  future  with  higher  grid  resolution  to  obtain  better  results.  Also,  the 
dispersion  model  simulations  will  be  extended  to  cover  the  entire  period  of  the  field 
experiment. 


4.  Trajectory  model 

4.1  Model  description 


This  trajectory  model  uses  interpolated  three-dimensional  wind  to  calculate  the 
displacement  of  an  air  parcel  during  each  time  step.  Height,  given  as  pressure  coori^ate, 
is  converted  to  meters  above  the  sea  level  by  using  the  hydrostatic  equation.  Either 
forward  trajectory  or  backward  trajectory  can  be  calculated  with  this  progr^.  For  large- 
scale  modeling,  for  which  a  rectangular  grid  is  no  longer  applicable  because  of  a 
projection  problem,  the  model  uses  latitudinal  and  longitudinal  to  identify  gnd  locations. 
The  model  also  considers  the  effect  of  terrain  on  the  parcel  trajectory  calculations. 


Meteorological  data  for  the  research  was  generated  by  the  Fifth-Generation  Penn 
State/NCAR  Mesoscale  Model  (MM5).  Wind  components  of  u,  v,  and  w,  pressure,  and 
temperature  on  all  sigma  layers  at  each  grid  point  were  used.  Terrain  information  such  as 
latitude,  longitude,  and  height  are  given  at  each  grid  point.  Grid  interval  is  set  to  31^. 
Lambert  conformal  projection  was  used  for  domain  setting.  The  projection  of  the 
spherical  earth  onto  a  plane  causes  each  grid  points  to  have  different  latitudinal  an 
longitudinal  spacing  despite  a  constant  distance  spacing  in  km.  The  trajectory  model 
interpolates  from  four  vertices  of  a  grid  cell. 

A  flow  chart  of  the  three-dimensional  trajectory  program  for  the  meteorological  data  is 
given  in  Figure  14.  Considering  forward  trajectory,  the  program  reads  the  initial  location 
(x,  y  p)  of  a  trajectory.  Meteorological  data  such  as  wind  components,  pressure,  and 
temperature  are  read  from  two  data  sets.  Grid  information  such  as  latitude,  longitude,  and 
teirain  altitude  are  also  given.  The  proportional  spatial  location  (xd,  yd  pd)  of  a 
trajectory  with  respect  to  the  south  west  comer  of  a  grid  cell  is  determined  from  the 
subroutine  LOCATE.  The  proportional  spatial  locations  have  values  between  0  and  1. 
Using  the  grid  wind  data  of  time  t,  xd,  yd,  and  pd,  the  wind  at  the  current  ^^tory 
location  can  be  calculated.  The  interpolation  is  performed  at  the  subroutine  INTRPL.  to 
the  subroutine  DISPL,  the  displacement  of  the  air  parcel  is  calculated  using  the 
relationship:  distance  (m)  =  time  (s)  *  velocity  (m/s).  To  calculate  backward  trajectoiy, 
the  time  step  is  set  to  be  negative. 

It  is  possible  that  the  calculated  new,first  guess,  trajectory  location  (x2,  y2,  p2)  is  located 
under  the  terrain  surface.  If  p2  is  larger  than  the  surface  pressure  at  location  (x2,  y2),  the 
trajectory  is  considered  to  be  under  the  surface  and  two  treatment  options  are  then 
considered  First,  the  surface-crossing  problem  is  treated  by  assuming  that  there  is  no 
vertical  wind  at  the  surface.  Therefore  the  trajectory,  once  intercepts  the  surface,  it  moves 
along  the  terrain  surface  with  only  horizontal  wind.  The  other  treatment  assumes  that  the 
air  parcel  is  reflected  by  the  surface  at  the  point  where  it  meets  with  the  surface. 
Subroutine  REFLECT  handles  the  reflection  process. 


Subroutine  LOCATE  gives  the  proportional  spatial  location  (xd,  yd,  pd)  at  the  new,  first 
guess,  trajectory  position.  With  wind  (u2,  v2,  w2)  and  spatial  location  (xd,  yd,  pd),  the 
interpolated  wind  (up2,  vp2,  wp2)  can  be  calculated  at  this  trajectory  position.  The 
average  of  (upl,  vpl,  wpl)  and  (up2,  vp2,  wp2)  gives  the  averaged  interpolated  wind  and 


Figure  14. 
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I 
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- 
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it  represents  the  wind  at  the  given  time  step  of  the  trajectory  (urn,  vm  wm).  Through  the 
subroutine  DISPL,  the  real  next  step  trajectory  position  (x2,  y2,  p2)  can  be  calculated 
using  the  averaged  interpolated  wind  (urn,  vm,  wm).  If  p2  is  greater  that  the  surface 
presiire  at  (x2,  y2),  reflection  will  be  carried  out  if  the  reflection  index  is  on.  Otherwse 
the  air  parcel  keeps  following  the  terrain  surface  after  it  meets  the  surface.  This  whole 
process  repeats  for  each  time  step  until  the  parcel  is  out  of  the  domain  or  the  time  expires. 

Figure  15  is  an  example  of  a  trajectory  simulated  using  the  trajectory  model.  The  initial 
time  is  07Z  10  August  1992.  The  result  is  for  the  three-dimensional  trajectory  with  a  one- 
minute  time  step,  using  original  wind  data  from  the  MM5  data  and  considering  reflection 
when  the  trajectory  meets  the  surface.  The  three  different  trajectones  m  the  top  plot 
reflect  the  three  different  initial  heights.  These  initial  heights  are  in  the  range  of  the 
effective  plume  height  from  the  source  during  the  experiment  day.  The  contours  ^ow 
surface  height  from  sea  level.  The  lower  half  of  Figure  15  shows  the  vertical 
displacement  along  the  three  trajectories.  The  dashed  lines  represent  the  trajertones  an 
the  solid  lines  represent  the  surface  height  along  each  trajectory.  From  this  figure,  not 
only  horizontal  variations,  but  also  vertical  variations  are  seen. 


4.2.  Sensitivity  tests 

The  parameters  that  affect  the  trajectory  calculations  are  mostly  related  to  wind  direction 
and  speed.  Interpolation  of  the  time  interval  and  imtial  height  are  important  parameters 
when  deciding  the  pattern  of  trajectories.  Determining  whether  the  considered  trajectory 
is  two-  or  three-dimensional  is  also  an  important  factor.  For  instance,  choosing  different 
surfaces  such  as  isentropic  surface,  iso-sigma  surface,  or  three-dimensional  trajectory  will 
result  in  different  trajectories.  The  present  study  compared  the  results  between 
trajectories  following  an  iso-sigma  surface  and  trajectones  considering  vertical  wind  (3- 
dimensional).  The  sensitivity  tests  are  performed  by  changing  one  of  the  following 
factors  For  the  wind  speed  factors,  the  model  was  run  with  three  different  speed 
categories,  i.e.,  the  original,  +10%,  and  -10%  of  the  modeled  speed  magmtudes, 
respectively.  For  example,  with  the  original  vertical  wind  speed,  three  trajectones  were 
calculated  by  only  changing  the  horizontal  wind  speed  among  onginal,  +10,  and  -lO/o  ot 
the  modeled  speed.  The  same  procedure  was  repeated  with  +10%  vertical  wind  speed  and 
-10%  vertical  wind  speed.  The  sensitivity  test  includes  all  factors  such  as  dimension, 
integration  time  steps  and  initial  heights.  We  can  then  figure  out  the  sensitivity  of  the 
trajectory  model  (Table  2). 

Table  3  shows  the  results  of  the  trajectory  ensemble  with  the  same  starting  time  but 
different  parameters.  The  trajectories  starting  at  07Z  August  1992  are  generally  toward 
the  northeast  but  we  can  see  the  trajectories  have  different  paths  depending  on  the  c^e. 
This  implies  that  changing  parameters  affects  the  trajectory  to  limited  extent.  For  other 
initial  times,  similar  results  were  obtained. 

To  estimate  the  results  of  a  trajectory  model,  the  present  study  uses  tracer  potential, 
which  has  been  developed  by  Koracin  et.  al.  (1999).  The  tracer  potential  is  given  by: 


Pr««8ur®(hPo)  Pr6«»ura(hPa)  Pro«8ur»(hPa)  LAXmJDE 


TRAJECTORIES  WITH  DIFFERE^fT  INITIAL  H0GHTS 


Figure  15.  Trajectories  of  air  parcels  released  at  three  different  heights  (925  mb,  930  mb, 

935  mb)  at  07Z 10  August  1992  from  the  source. 


Table  2  Parameters  considered  in  the  sensitivity  tests 


Parameters  for  sensitivity  test 

Caises 

Dimension  of  trajectory 

iso-sigma  trajectory 

Three-dimension  trajectory 

Time  step 

1  minute 

10  minute 

30  minute 

Initial  heights 

935mb 

930mb 

925mb 

Horizontal  wind 

Original  wind 

+10%  from  original  horizontal  wind 
-10%  from  original  horizontal  wind 

Vertical  wind 

Original  wind 

+10%  from  original  vertical  wind 
-10%  from  original  vertical  wind 

Tables 


Parameter  numbers,  parameters,  and  corresponding  tracer  potential  values 


Parameter  Number 

Parameter* 

Tracer  Potential 

1 

fQ001D2R0H0V0 

0.609 

2 

10001 D2R0H1  VO 

0.665 

3 

ID0O1D2ROH2VO 

0.841 

4 

f0101D2R0H0V0 

0.538 

5 

f0101D2ROH1VO 

0.648 

6 

10101 D2R0H2V0 

0.669 

7 

10201 D2R0H0V0 

0.682 

8 

10201 D2R0H1  VO 

0.770 

9 

10201 D2R0H2V0 

0.848 

10 

10001 D3R0H0V0 

1.305 

11 

f0001D3ROH0V1 

1.288 

12 

10001D3ROH0V2 

1.456 

13 

10001 D3R0H1  VO 

1.335 

14 

10001 D3R0H1 VI 

1.188 

15 

10001 D3R0H1V2 

1.401 

16 

10001 D3R0H2V0 

1.202 

17 

10CX)1D3R0H2V1 

1.743 

18 

10001 D3R0H2V2 

1.361 

19 

fD101D3R0H0V0 

0.986 

20 

fQ101D3ROHOV1 

1.427 

21 

10101 D3R0H0V2 

1.443 

22 

101O1D3ROH1VO 

1.283 

23 

1D101D3R0H1V1 

1.350 

24 

10101 D3R0H1V2 

1.337 

25 

1D101D3R0H2V0 

1.145 

26 

10101 D3R0H2V1 

1.319 

27 

10101 D3R0H2V2 

1.539 

28 

10201 D3R0H0V0 

1.303 

29 

10201 D3R0H0V1 

1.470 

30 

f0201D3R0H0V2 

1.500 

31 

t0201D3R0H1V0 

1.370 

32 

10201 D3R0H1 VI 

1.649 

33 

10201 D3R0H1V2 

1.595 

34 

f0201D3R0H2V0 

1.477 

35 

f0201D3R0H2V1 

1.593 

36 

10201 D3R0H2V2 

1.571 

*  The  pflTamptPT  ‘f  means  the  considering  trajectories  are  forward  trajectories.  The  following  two  digits 
represent  initial  hei^t  as  935mb,  930mb,  and  925mb  for  00,  01,  and  02,  respectively.  The  next  two  digits 
representing  time  step  ,in  minutes,  e.g.  ‘01’  means  Imin.  The  number  next  to  the  character,  ‘R’  represents 
either  the  trajectory  considers  reflection  at  flie  surface  (‘  1  ’)  or  not  (‘0’).  In  the  present  research,  no  reflected 
trajectories  are  considered  and  therefore  the  number  is  ‘O’.  Digits  behind  the  characters  H  and  V  are 
about  wind  changes.  ‘0’  reixesents  interpolated  wind  using  original  data,  ‘1’  represents  wind  q)eed  is  10% 
more  than  the  original  interpolated  value,  ‘2’  represents  wind  is  10%  less  than  the  original  data. 


where  Ci  is  the  concentration  at  receptor  I,  and  Ri  is  the  distance  between  receptor  and 
trajectory.  A  is  a  constant,  which  provides  a  convergent  solution  when  the  trajectory 
passes  over  the  receptor.  The  unit  for  tracer  potential  is  (femtoliter/liter/km). 

Koracin  et  al.  pointed  out  that  most  model  evaluation  schemes  apply  a  dispersion  model 
to  the  wind  fields  and  compare  dispersion  estimates  with  measured  concentration, 
assuming  that  measured  concentration  is  real.  However,  since  the  wind  field  simulated  by 
models  over  complex  terrain  are  usually  of  poor  quality  it  is  important  to  evaluate  the 
wind  fields  for  dispersion  modeling. 


In  calculating  the  tracer  potential  of  a  trajectory  position,  concentration  at  each  receptor 
was  weighted  by  both  the  magnitude  of  the  receptor’s  concentration  and  by  the  inverse  of 
the  distance  between  the  receptor  and  the  trajectory  position.  This  implies  that  if  we 
compare  the  weights  of  two  receptors  with  the  same  concentration,  a  receptor  closer  to  an 
object  trajectory  has  a  larger  weight  than  the  other  receptor.  This  also  implies  that  if  we 
compare  the  weights  of  two  receptors  with  same  distance  from  the  object  trajectory,  the 
receptor  having  a  larger  concentration  has  a  larger  weight  than  the  other  receptor. 


Tracer  potential  (TPo)  at  the  source  is  put  as  reference  value. 

n  n 

where  Q  is  concentration  at  receptor  i,  and  Rj  is  distance  between  the  source  and 
receptors. 


A  tracer  potential  larger  than  the  reference  tracer  potential  implies  success  of  wind  field. 
By  comparing  tracer  potential  with  the  reference  tracer  potential,  we  can  quantify  the 
success  of  the  wind  field  itself  without  relative  comparison  to  other  wind  fields.  Another 
advantage  of  the  tracer  potential  method  is  that  even  if  sampling  sights  are  out  of  the 
domain  considered  for  numerical  modeling,  values  fi-om  these  sampling  sites  can  be  used. 
In  other  words,  there  is  no  waste  of  observational  data  fi’om  field  program. 


The  results  of  sensitivity  tests  by  changing  parameters  are  given  in  Table  3.  The  results 
are  shown  in  Figure  16.  The  sensitivity  test  shows  that  three-dimensional  trajectories 
represent  better  wind  fields  better  than  two-dimensional  trajectories.  For  the  wind  data, 
using  the  original  MM5  data  gives  almost  the  same  values  at  the  initial  level  of  93 5 mb 
and  925mb.  But  at  the  initial  level  of  930mb,  the  tracer  potential  value  is  26%  lower  than 
other  levels,  which  implies  that  the  contribution  of  930mb  is  less  than  other  levels. 


For  the  three-dimensional  cases,  the  925mb  initial  height  gives  the  highest  tracer 
potential  value  as  1.503([fl/l]*[hr/km])  and  the  lowest  standard  deviation  value  of  0.113. 
For  the  two  lower  levels,  the  averages  and  standards  deviations  of  tracer  potentials  for 
three-dimensional  trajectories  don’t  have  significant  differences.  The  average  value  of 
tracer  potential  values  including  all  parameters  is  1.394([fl/l)*[hr/km]),  and  the  standard 


Positive  tracer  potential  ([fl/l]*[hr/kml) 


Sensitivity  Tests  by  Tracer  Potential  (Imin  time  step) 


Figure  16.  Sensitivity  tests  of  the  trajectory  model  by  using  tracer  potential 


Table  4.  Average  and  Standard  deviation  of  tracer  potential 


Average 

STD  (AlO 

STD  (3D-only) 

935mb 

1.372 

0.334 

0.197 

930mb 

1.355 

0.346 

0.165 

925mb 

1.503 

0.349 

0.113 

Total 

1.394 

0.341 

0.116 

deviation  is  0.166.  It  seems  that  the  effective  plume  height  of  the  tracer  in  925mb 
contributes  more  than  at  other  levels. 

4.3  Conclusions 

It  is  shown  that  a  three-dimensional  wind  field  is  necessary  in  trajectory  analysis. 
Reflection  at  the  surface  has  no  major  effect  on  results.  At  least  in  this  case,  the  accuracy 
of  wind  speed  does  not  contribute  greatly  to  the  tracer  potential  analysis.  However, 
correct  wind  field  data  is  always  beneficial. 

This  trajectory  has  many  applications.  It  will  be  used  to  determine  the  aerosol  origins  in 
various  parts  of  Texas  to  see  if  the  chemical  content  of  aerosol  has  any  relationship  with 
rain-cloud  formation.  We  are  also  planning  to  use  this  model  to  trace  the  origin  of  air  in 
the  Antarctic  during  years  with  El  Nino  or  years  with  La  Nina  to  see  if  there  are  any 
significant  differences  between  them.  This  may  help  in  analyzing  ancient  climate 
changes. 

5.  Summary 

Three  scientific  papers  and  a  master’s  thesis  are  under  preparation  from  this  project.  A 
higher-order  closure  model  suitable  for  dispersion  simulation  in  complex  terrain  has  been 
developed.  The  model-simulated  features  in  CBL-dispersion  with  different  source 
heights  have  been  performed.  The  results  are  in  very  good  agreement  with  those  from  the 
tank  experiments.  It  is  shown  that  the  advection  terms  in  the  equations  for  the  second 
order  moments  that  involve  concentration  are  important  in  simulating  the  essential 
features  of  dispersion.  The  model  shows  good  agreement  when  applied  to  the  MOHAVE 
1992  field  data. 

Comparison  of  wind  fields  predicted  by  both  the  MM5  and  the  MIUU  mesoscale  models 
does  not  give  conclusive  statistical  answers  on  which  model  performs  better.  It  is 
probably  better  to  use  a  trajectory  model  to  evaluate  the  wind  field  predicted  by  various 
models  (Koracin  et  al.,  2000). 


APPENDIX 


A.  THE  MIUD  MODEL 

The  MIUU  model  is  a  three-dimensional,  hydrostatic,  incompressible,  higher  order 
turbulence  closure  model.  The  model  was  developed  at  the  Department  of  Meteorology, 
Uppsala  University  (MIUU)  Sweden.  A  terrain-influenced  co-ordinate  system,  Ti-system, 
is  used  in  the  model.  The  basic  equations  of  the  model,  transformed  to  this  ti  co-ordinate 
system  are  presented  below. 
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where  W  is  the  mean  vertical  wind  speed  in  the  terrain  following  coordinate  system,  ti. 
For  more  information  see  e.g.  Enger  (1986),  Tjernstrdm,  (1988),  Andren  (1990) 


B.  LENGTH  SCALE  FORMULATION 

The  turbulent  length  scales  used  in  the  model  scales  are  parameterized  “J  a  similar 
fashion  with  Enger  (1986,  1990)  and  included  the  modification  by  Andren  (1990).  It  is 
supposed  that  all  the  length  scales  are  to  be  proportional  to  each  other.  This  means  one 
only  needs  to  determine  one  length  scale  (X),  which  is  sufficiently  general  to  handle 
various  types  of  flow  in  the  atmosphere. 
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A  smooth  transition  is  obtained  between  (Bl)  and  (B2)  by  taking  the  minimum  of  the  two 
expressions. 


Neutral  and  stable  boundary  layer; 
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At  the  top  of  the  boundary  layer  these  formulation  are  limited  by  using  the  minimum  of 
the  stability-dependent  formulation. 

®B 

where  m  is  the  Brunt-Vaisala  frequency  and  the  neutral  form  of  equation  (B3).  For  the 
constant  Cb  a  value  of  0.5  is  used.  For  constant  Bi  a  value  of  21.33  is  used,  which 
incorporate  the  parameterization  of  ‘wall  effect’  at  surface  as  introduced  by  Andren 
(1990).  For  more  information  on  the  calculations  of  the  constants  see  Enger  (1986, 
1990a)  and  Andren  (1990). 
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NOMENCLATURE 


Symbols  that  are  defined  in  the  text  are  not  included  in  this  list 


C 

c' 

D 

f 

gi 

K 

L 

Prt 


mean  concentration 
fluctuating  part  of  concentration 
molecular  destruction  term 
Coriolis  parameter 
acceleration  of  gravity,  gi  =  (0,0,-g) 

general  notation  for  scalar  eddy  transfer  coefficient,  e.g.  Kr,  and  Km  are  eddy 
transfer  coefficients  for  heat  and  momentum,  respectively 

Monin-Obukhov’s  length 

Prandtl  number,  P^  =  w’0'^ 


p'  fluctuating  part  of  the  pressure 

mean  double  turbulent  energy 
Q  source  strength 

R  mean  specific  humidity 

Sc  Sinks  and  sources 

s  Model  top 

U,V,W  mean  wind  speed  in  x,  y,  and  ti-  direction,  respectively 
Ug,Vg  geostrophic  wind  in  x-  and  y-direction,  respectively 
Ui  mean  wind  in  tensor  form 

u'yy  fluctuation  part  of  the  wind  speed  in  x,  y,  and  z-  direction,  respectively 
u’i  fluctuation  part  of  the  wind  speed  in  tensor  form 


friction  velocity 


u. 

w .  convective  velocity  scale,  =  (z, 

x, y,z  coordinate  in  west-east,  south-north,  and  vertical  direction,  respectively 

Xj  coordinate  in  tensor  form 

Zg  terrain  height 

Zj  convective  boundary  layer  height  defined  as  the  height  of  the  minimum  heat 

flux 

ai,  a2,pi,  pH  constants 
X  characteristic  turbulent  length  scale 

Z  —  Zg 

Ti  terrain  influenced  coordinate  system  (=  5—  ) 

p  density 

Or  radiative  heating/cooling  rate 

0  mean  potential  temperature 

0’  fluctuating  part  of  the  potential  temperature 

n  a  scaled  pressure,  the  so-called  Exner  function 


